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ABSTRACT 

The recent detection of a giant exosphere surrounding the warm Neptune GJ 436 b has shed new light on the evapora¬ 
tion of close-in planets, revealing that moderately irradiated, low-mass exoplanets could make exceptional targets for 
studying this mechanism and its impact on the exoplanet population. Three HST/STIS observations were performed in 
the Lyman-a line of GJ 436 at different epochs, showing repeatable transits with large depths and extended durations. 
Here, we study the role played by stellar radiation pressure on the structure of the exosphere and its transmission 
spectrum. We found that the neutral hydrogen atoms in the exosphere of GJ 436 b are not swept away by radiation 
pressure as shown to be the case for evaporating hot Jupiters. Instead, the low radiation pressure from the M-dwarf 
host star only brakes the gravitational fall of the escaping hydrogen toward the star and allows its dispersion within 
a large volume around the planet, yielding radial velocities up to about -120kms“^ that match the observations. We 
performed numerical simulations with the EVaporating Exoplanets code (EVE) to study the inffuence of the escape 
rate, the planetary wind velocity, and the stellar photoionization. While these parameters are instrumental in shaping 
the exosphere and yield simulation results in general agreement with the observations, the spectra observed at the 
different epochs show specific, time-variable features that require additional physics. 
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1. Introduction 

1.1. Atmospheric escape 

More than 40% of the known exoplanets orbit extremely 
close to their star (^0.1 an), raising many questions about 
the state of their atmosphere and its interactions with the 
host star. Theoretical studies indicate that the intense stel¬ 
lar X-ray and extreme ultraviolet energy input into the 
lower thermosphere of hydrogen-rich planets can lead to a 
hydrodynamic expansion of the upper gas layers. If the flow 
of expanding hydrogen reaches high velocities or high alti¬ 
tudes, it can escape the gravitational pull of the planet in 
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Transit spectroscopy in the ultraviolet has proved an in¬ 

valuable tool for studying the outer regions of an evaporat¬ 
ing planet. Because of its expansion, the atmosphere pro¬ 
duces a much deeper absorption than the planetary disk 
when observed in the lines of elements that are abundant 
at high altitudes and/or associated with strong electronic 
transitions. Observations in the bright stellar Lyman-a line 
of neutral hydrogen have led to the detec tions of evap¬ 
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et al.|2Q12 ). Heavier metals and ions can be carried upward 
through collisions with the expanding flow of hydrogen, 
which was confirmed by detection of several s uch species at 
high altitudes around HD 20945 8 b G+, Vidal-Madiar 
et al. 2004 Linsky et al. 2010| Mg^, | Vidal-Madiar et al. 


20131 and tentat ively Linsky et al.||201Q Ballester fc 
Ben-Jaffel 2015|) and HD 18 9733 b (O^ and possibly G+, 


Ben-Jafi'el & Ballester 2013), and allowed detection of an 
extended exosphere around the hot-Jupiter WASP-12 b 
through transit observations in the ne ar-UV (Mg+, Ee+; 


Eossati et al.||2010 , Haswell et al.||2012 ). 

Observations at high-resolution in the UV have also been 
used to probe the structure of these extended exospheres, 
revealing that they are shaped by interactions with the host 
star, such as radi ation pressure, stellar wind, and magnetic 
interactions (e.g. . Llama et al.||2011[_Bourrier fc Lecavelier 
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Observations at 


HD 189733 b led 
further to the first detection of temporal variations in the 
physical properties of an extended exosphere, which prob¬ 
ably come from the energetic and/or dynamical influence 
of the host star ( Lecavelier des Etangs et al.|2012 Bourrier 
et al.|[2QT3 ). 


Many close-in exoplanets are expected to show signs of 
evaporation. The impact of this phenomenon on the forma- 
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moderate escape rate, leading to the loss of a few per¬ 
centage points of the planet mass over the lifetime of the 
system, low-density planets like mini-Neptunes or super- 
Earths with a large volatile envelope may be significantly 
eroded by evaporation. This scenario is supported by trends 
in the e xoplanet population: the lack of sub-Jupiter size 


planets dLecavelier des Etang! 


2009 Beauge &: Nesv orny |2013) and the l arge number of 
close-in super-earths ( [Howard et ah [2012 ) are consistent 


2007 Davis & Wheatley] 


with the formation of rocky planetary remnants through the 
evaporation of their hydrogen/helium envelope. It is thus 
crucial to observe the atmosphere of close-in Neptune-mass 
planets to understand in which conditions their evapora¬ 
tion can lead to the formation of atmosphereless remnant 
cores like CoRoT-7b and Kepler-10 b. 


1.2. The case of GJ436 b 


In this framework, the warm Neptune GJ436 b (Butler 
|et al.|[2QQ4{ |Gillon et al.||2007| ) is a good target for transit 
observations in the Lyman-a line, as it orbits a bright star 
(V=10.7) whose close proximity to the Earth (d=10.3pc) 
limits the effects of interstellar medium absorption in the 
core of this line ( Ehrenreich et alT]|2011 ). The low mass of 
GJ436 b (Mp=23.1 MEarth)^ the m oderate irradiation from 
its M-dwarf host star (a=0.029 an; |Ehrenreich et ar]|2015|), 
and the pu zzling eccentricity of its orbit (e=0.16, Beust 


et al. 2012 and references therein) offer the oppo r tunity 
to study a new evaporation regime. |Kulow et al. (2014) 


reported the detection of a deep absorption signature in the 
stellar Lyman-a line starting near the end of the optical 
transit, thought to be caused by the occultation from 
a large cometary tail of neutral hydrogen H^. However, 
interpretation of these observations was misled by the use 
of an inaccurate transit ephemeris and by the lack of an 
out-o f-transit reference for th e flux in the stellar Lyman-a 
|2015[). U s ing tw o additional HST 
T et c 


line ( Ehrenreich et al. 
observations, renreicJ 


2. HST/STIS observations in the Lyman-a line 

We used the three existing transit data sets of GJ 436b ob¬ 
servations in the Hi Lyman-a line (1215.6702 A), all taken 
with similar settings with the Space Telescope Imaging 
Spectrograph (STIS) instrument onboard the Hubble Space 
Telescope (HST). The HST data consist of time-tagged 
spectra obtained with the G140M grating, with a spectral 
resolution of ^ 20 kms“^ at 1215.67 A. Data was reduced 
as described in [Ehrenreich et a'r]|2Q15[ Throughout the pa¬ 
per we use the term transit to refer to the transit of the 
extended exosphere of GJ 436 b, while we explicitely speak 
of optical transit to refer to the occultation caused by the 
planetary disk alone. Each data set was obtained through 
four HST orbits performed at different phases near the op¬ 
tical transit. During Visit 1 (December 2012 ; [Kulow et al. 
2014), one exposure was taken before the optical transit, 
one during and two after. V isits 2 and 3 (June 2013 and 


2014; Ehrenreich et al.|2015 ) have two exposures before the 
optical transit, one during and one after. The tim es of mid- 
exposures were calculated with the ephemeris in Knutson 
et al.|[2QTTj and are reported in Table 


3. Reconstruction of the Lyman-a line 

3.1. Method 

The observed Lyman-a line is affected by absorption from 
the interstellar medium (ISM), which needs to be corrected 
for in order to estimate the intrinsic stellar Lyman-a line 
profile as seen from the planet. It is used to calculate the 
radiation pressure on the escaping exosphere (which is de¬ 
pendent on the stellar flux), and as a reference to calcu¬ 
late the theoretical transmissi on sp ectra of the exosphere 
as seen with HST/STIS (Sect. |4.3| ). The intrinsic line can 
be reconstructed directly using the STIS observations, and 
we applied the same method to all visits independently 


_ ^ (2015) revealed a much 

deeper Hi absorption that starts several hours before the 
optical transit (^2 - 4h) and that may remain detectable 
for more than 20 hours afterward. GJ 436 b is thus not only 
trailed by a cometary tail, but also surrounded by a giant 
coma of neutral hydrogen extending tens of planetary radii 
in front of the planet. 

Our goal is to investigate the effects of radiation pressure 
and the basic physical parameters of an evaporating 
exoplanet on the structure of GJ 436 b hydrogen exo¬ 
sphere and its transmission spectrum. To this aim we 
reconstructed the intrinsic Lyman-a line of the host star 
and compared the results of 3D numerical simulations 
of extended atmospheres based on the EVaporating 
Exoplanets code (EVE) with the three sets of transit 
observations. Observations are presented in Sect. and 
the reconstruction of the intrinsic stellar Lyman-a line is 
explained in Sect. The EVE code is described in Sect.|^ 
and used in Sect. ^ to study the properties of GJ 436 b 
exosphere in the low radiation pressure field of its host 
star. In Sect. we discuss the part played by radiation 
pressure in the observed features and investigate several 
mechanisms that could also intervene. We summarize our 
conclusions in Sect. 0 


(Wood et al. 2005 [Ehren reich et al. 2011, Erance et al. 
20121 [Bourrier et al. |2013p . We calculate a model profile of 


the stellar Lyman-a line, which is absorbed by the interstel 
lar hydrogen and deuterium, convolved by STIS line spread 
function (LSE), and then compared with an observed ref¬ 
erence spectrum. Eor each visit we used as reference the 
spectrum obtained during the first HST orbit, assuming it 
is unabsorbed by the exosphere of GJ 436 b. We used the 
Bayesian Information Criterion (BIC; 1^ Wit et al.[[M^ 
to prevent over-fitting in our search for the best moael for 
the Lyman-a and LSE profiles. Trying out different models 
for the stellar line and the LSE, we found that in all visits 
the lowest BICs were obtained with a single Voigt profile 
for the stellar emission line, and a single Gaussian profile 
for the LSE. The free parameters of the fit are thus the 
heliocentric velocity of the star 7 *, the Lyman-a line pro¬ 
file parameters (maximum stellar flux at 1 an fpea/c(lau), 
Doppler width Ayjj and damping parameter a damp) ^ the 
width of the LSE profile (Tlsf^ and the ISM parameters 
for neutral hydrogen (column density logio ^(Hl), Doppler 
broadening parameter 5(Hl), and velocity relative to the 
star 7 (Hi)/^). Although there may be multiple ISM com¬ 
ponents along the line of sight (LOS) toward GJ 436, we 
assumed their absorption of the stellar Lyman-a profile can 
be interpreted as that of a single component with average 
properties. We thus modeled the ISM opacity as the com¬ 
bination of two Voigt profiles for the atomic hydrogen and 
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Time 

Pre-transit 

Transit 

Post-transit 

Visit 1 

- 

[-01:55 ; -01:30] 

[-00:43 ; 00:05] 

[00:52 ; 01:41] [02:28 ; 03:17] 

Visit 2 

[-03:23 ; -02:55] 

[-02:01 ; -01:27] 

[-00:26 ; 00:09] 

[01:10 ; 01:45] 

Visit 3 

[-03:26 ; -02:58] 

[-02:00 ; -01:25] 

[-00:24 ; 00:10] 

[01:11 ; 01:46] 


Table 1. Log of GJ 436b transit observations. Time is given in hours and minutes, and counted from the center of the 
optical transit. 


deuterium, separated by about 0.33 A. We used a Dl/Hl 
ratio of 1.5x10“^ Hebrard fc Moos|2QQ3 Linsky et al. 
|2QQ6[ ) and assumed there is no turbul ent broadening with 
6(Hl) = b(D l)/\/2 jWood et al.||2QQ5 ). 

A preliminary analysis of variations was used to lo¬ 
cate the best model, and revealed no multiple local min¬ 
ima in the parameter space. We then refined the best-fit 
parameter values and calculated their uncertainties using 
a Metropolis-Hasting Markov chain Monte Carlo (MCMC) 
algorithm. To better sample the posterior distribution in 
the case of non-linear correlations between parameters, we 
applied an adaptive principal component analysis to the 
chains and jumpe d the parameters in an uncorrelated space 
( Diaz et ^|2Q14 ). The step size was adjusted to ensure an 
acceptation rate of about 25%, and the system was ana¬ 
lyzed with chains of 5x10^ accepted steps. The comparison 
between the model spectra and the observed spectra was 
performed over the entire spectral line, excluding the core 
where ISM absorption is saturated and noise from the air- 
glow correction is high (the spectral range of the fits are 
given in Table [^. The medians of the posterior probability 
distributions were chosen as the final best-fit values for the 
model parameters, and their la uncertainties were obtained 
by finding the intervals on both sides of the median that 
contain 34.15% of the accepted steps. 


3.2. Results 

The best-fit reconstructed Lyman-a stellar profiles are dis¬ 
played in Fig. with the corresponding model parame¬ 
ters given in Table The reconstructed Lyman-a lines for 
Visits 2 and 3 are remarkably similar, with all fit parameters 
consistent within their la uncertainties. On the other hand 
we found spurious values for the best-fit parameters in Visit 
1, with the ISM absorption profile anomalously deep and 
strongly blue-shifted. This is consistent with a lower fiux 
in the blue wing of the Lyman-a line caused by absorp¬ 
tion from hydrogen escaping the planet, and confirms that 
the spectrum obtained during the first exposure in Visit 
1 is already absorbed by the giant exosphere of GJ 436 
b (Ehrenreich et al. 2015). Without a good knowledge of 
both wings of the line, it is not possible to reconstruct the 
intrinsic stellar line for Visit 1. However the atmospheric 
absorption occurs in a localized spectral rang e in the blue 
wing of the line, and |Ehrenreich et alT]|2Q15| showed that 
the Lyman-Of flux measured at larger velocity blue shifts, 
and over the whole red wing of the line, is very stable be¬ 
tween the three visits, which is in agreement with our line 
reconstructions in Visits 2 and 3. Hereafter we thus use the 
Lyman-Of line profile obtained for Visit 2 as a proxy for 
Visit 1. 

The best-fit parameters we obtained for Visits 2 and 3 cor¬ 
respond to heliocentric radial velocity for the ISM hydro¬ 


gen of -1.7T2.1 kms“^ and -l.QT2.7k ms~^, respectively. 
Using the LISM Kinematic Calculatoip (Redfield & Linsky 
2008), we identified three ISM clouds near (less than 
20 g fhe line-of-sight toward GJ 436 with radial velocities 
of 0.37±1.38kms“^ (LIC cloud), 2.33±0.76kms“^ (Leo 
cloud), and 5.81±0.74kms“^ (NGP cloud). The ISM ab¬ 
sorption derived from our reconstruction of the Lyman-a 
line thus corresponds more likely to the LIC cloud (<lcr), 
and may be consistent with the Leo cloud (<2cr). 


4. EVE: the EVaporating Exoplanets code 

EVE is a 3D numerical code that we developed to calculate 
the structure of an exoplanet upper atmosphere and the 
transmission spectrum of the species it contains. The di¬ 
rect comparison of these spectra with observational data al¬ 
lows the inferrance of several physical properties of both the 
planetary atmosphere and the star, and the code has been 
used in this way to study the neutral hydrogen e xospheres 


of the hot Jupiters HD 2Q9458 b and HD 189733b (Bourrier 


and the super-Earth 55 Cnc 


(Ehrenreich et al. 

2012 

), as 

. magnesium 

popu. 

lations in 


the th ermosphere/exosphere of H D 2Q9458b (Bourrier et al. 
2014). The EVE code was used in Ehrenreich et al.l(|2015|) 
to perform a preliminary fit of the three combined ooserva- 
tions of GJ 436b neutral hydrogen exosphere. In the present 
paper, we further investigate the effects of stellar radiation 
pressure, stellar photo-ionization, and the planetary escape 
properties on the structure of the exosphere and its trans¬ 
mission spectrum. The following subsections offer a syn¬ 
thetic description of the model used in this frame. 

4.1. Bottom layers of the atmosphere 

The upper atmosphere is divided in two different regimes, 
joined at a transition altitude set at i?trans=4 J^pianet, 
approximately the mean altitude of the Roche lobe. We 
assume that the atmosphere within the Roche lobe, which 
is by definition under the gravitational influence of the 
planet, keeps a global cohesion and can be described 
analytically with a radial density profile corresponding to 
a hydrostatic equilibrium. In the numerical simulation, 
the density profile of neutral hydrogen in this part of the 
atmosphere is reassessed every time step by adjusting its 
value at i^trans to match the hydrogen de nsity obtained in 
the upper atmospheric layers (see Sect. 4.2). Because of 
the eccentricity of the planetary orbit (e=0.16), the density 
profile thus varies over time. The gas filling the Roche 
Lobe mainly contributes to the absorption in the core 
of the Lyman-a line, which cannot be observed because 


http://sredfield.web.wesleyan.edu/ 
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Table 2. Parameters derived from the reconstruction of the Lyman-a line profile in Visits 2 and 3. The line could not 
be reconstructed in Visit 1 due to the lack of reference spectrum. 



Parameter 

Visit 2 

Visit 3 

Unit 

Fit range 

Blue wing 

[1214.64 - 1215.53] 

[1214.68 - 1215.53] 

A 


Red wing 

[1215.75 - 1216.70] 

[1215.75 - 1216.54]t 

A 

ISM (Hi) 

logic iV(Hi) 

18.01±0.07 

17.90io:i9 

cm“^ 


6(Hi) 



kms“^ 


7(Hi)/* 

-9.iii:l 

-9.2tl:l 

kms“^ 

LSF 

CTLSF 

0.97±0.14 

i.25i!^:if 

HST pixels 

Stellar line 

^peak (1 au) 

1.55-0.24 

1 ^9+0-56 

1.0Z_o.38 

erg cm“^ 



74±3 

69±5 

kms“^ 


ridamp 

0.09±0.02 

O.lOiaSi 



7* 

7.4±1.6 

8.2if:^ 

kms“^ 


EL 2 /a (Earth) 

2.05x10“^^ 

2.01x10"^® 

erg cm“^ 


EL^a(lau) 

0.90 

0.88 

erg cm“^ 


Note: F Lya is the flux integrated in the entire Lyman-o line and calculated at Earth distance and at 1 an from the star, 
f: An anomalous low-flux pixel at about 1206.OA was removed from the fit (see Fig.j^. 


of ISM absorption and geocoronal emission. Therefore, 
this part of the atmosphere has little influence on the 
observed spectra, except through the Lorentzian wings of 
the absorption p rofile (the so-called damping wings, e.g. 


Ben-Jaffel 2QQ8[ ). This is however excluded because no 
significant absorption signal was observed in the red wing 
of the Lyman-g line du ring the optical transit of GJ 436 b 


( Ehrenreich et al.]|2Q15 ) 


4.2. Upper layers of the atmosphere 

We use Monte-Carlo particle simulations to compute the 
dynamics of the escaping gas above i^trans- Neutral hydro¬ 
gen atoms are represented by meta-particles, and the total 
number of metaparticles launched every time step dt, set to 
5 mn, depends on the escape rate of neutral hydrogen along 
the eccentric orbit. To take the eccentricity into account, 
we assumed that GJ 436b is in the energy-limited regime of 
evaporation, and the escape rate is proportional to the stel¬ 
lar energy input into the atmosphere. We thus varied the 
escape rate as the inverse square of the distance to the star, 
setting the reference Mho at the semi-major axis. The parti¬ 
cles are released from the entire atmosphere and their initial 
velocity distribution, relatively to the planet, is the combi¬ 
nation of the upward bulk velocity of the planetary wind 
'^pi-wind and an additional thermal speed com ponent from 
a Maxwell-Boltzman veloci ty distribution (see Bourrier & 
Lecavelier des Etangs||2Q13 ). 

The dynamics of the particles is calculated in the stellar 
reference frame and results from the stellar and planetary 
gravities, the stellar radiation pressure, and the inertial 
force linked to the non-Galilean reference frame. Radiation 
pressure on an atom is the combination of a radiaj^ im¬ 
pulsion in the direction opposite to the star, caused by 
the absortion of a photon, followed by an isotropic impul- 


^ Hereafter, we use the term projected velocity to refer to the 
projection of a particle velocity on the star-Earth line of sight, 
and radial velocity to refer to its projection on the star-particle 
axis. Radiation pressure varies with the radial velocity of a parti¬ 
cle, while we measure spectra as a function of projected velocity. 


sion after the photon reemission. We process the cumulated 
radial impulses from all photons absorbed by a hydrogen 
atom during dt as a force proportional to stellar gravity 
( Lagrange et al.|[l998 ), and apply an equivalent number of 
velocity changes with constant speed 3.26 ms“^ (this value 
comes from the conservation of momentum) and random 
direction to reflect the isotropic impulsions. A crucial as¬ 
pect of radiation pressure is its dependence on a hydrogen 
atom radial velocity, as the number of photons it absorbs 
and reemits is proportional to the flux in the Lyman-a line 
Doppler-shifted at this velocity. The effect of radiation pres¬ 
sure thus strongly depends on the shape of the line profile, 
and the single peak profile obtained for GJ 436b (Fig.[^ im¬ 
plies that the strength of radiation pressure monotonously 
decreases from the core to the wings of the line, ie with 
increasing absolute radial velocities. 

The hydrogen atoms lifetime depends on the stellar photo¬ 
ionization rate Tion which varies as the inverse square of 
the distance to the star. Once ionized, hydrogen atoms are 
removed from the simulation. 

We do not take collisions into account in the upper atmo¬ 
spheric layers, as neutral hydrogen densities reach about 
4xl0^m“^ at i^trans foi* OU.T bcst-fit simulations. With our 
simulation settings this value corresponds to a Knudsen 
number (the ratio between the mean free path of the gas 
and the atmospheric scale height) of about 40, which is 
above the level of the exobase defined by a Knudsen of one, 
even if the neutral fraction of hydrogen in the upper atmo¬ 


sphere of GJ 436 b is about 0.1 ( [Knlow et al.|2014 
et al.|[Wl3a ). 


Koskinen 


4.3. Theoretical spectra 

Self-shielding for radiation pressure and photo-ionization is 
taken into account by calculating their respective opacities 
in the bottom and upper atmospheric layers. Similarly, 
EVE calculates the opacity along the star-Earth line of 
sight to obtain at each time step the theoretical absorption 
spectra in the Lyman-a line ( for details regarding the 
calculation of these opacities, see [Bourrier fc Lecavelier des 
Etangs 2013). During the optical transit, the stellar flux 
is absorbed at every wavelength by the planetary disk. In 
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Doppler velocity (km/s) 


Fig. 1. Lyman-a line profiles of GJ 436 for Visits 2 (top) 
and 3 (bottom). The yellow line shows the theoretical in¬ 
trinsic stellar emission line profile as seen by hydrogen 
atoms escaping the planetary atmosphere, scaled to the 
Earth distance. It also corresponds to the profile of the 
ratio P between radiation pressure and stellar gravity (val¬ 
ues reported on the right axis). The solid black line shows 
the Lyman-(a line profile after absorption by the interstel¬ 
lar hydrogen (1215.6A) and deuterium (1215.25A), whose 
cumulated profile is plotted as a dotted black line (ISM 
absorption in the range 0-1 has been scaled to the entire 
vertical axis range). The line profile convolved with the LSF 
(blue line) is compared to the observations shown as a red 
histogram. Above the horizontal dashed line (/5=1), radia¬ 
tion pressure would overcome stellar gravity and accelerate 
particles away from the star. However with (3 values below 
unity, radiation pressure from GJ 436 can only brake the 
gravitational deviation of particles toward the star. 


contrast, the neutral hydrogen atmosphere may contribute 
to the absorption on longer time scales and at specific 
wavelengths that depends on the dynamics of the gas and 
the effects of thermal and natural broadening. At a given 
time step in the simulation, the theoretical intrinsic stellar 
line profile is multiplied by the absorption profile of the 
planet and its atmosphere, calculated with a resolution 
/S.V =10kms“^ (more than twice the resolution of the 
STIS spectra). The resulting spectrum is then attenuated 
by the ISM absorption and conv olved with the LSF (using 
the fits obtained in Sect. 3.2), and interpolated to the 
observed velocities. To compare the simulation results with 
the STIS observations, we average the theoretical spectra 
calculated during the time window of each observation. 


5. Structure of the exosphere 

While mechanisms such as stellar wind or bow-shocks may 
play a role in shaping the exosphere of an exoplanet, their 
description requires strong assumptions and the introduc¬ 
tion of several free parameters in models of atmospheric es¬ 
cape. By contrast radiation pressure is known to be present 
with a well defined behavior, and we showed in Sect.j^that 
it is very well constrained from direct measurements of the 
stellar Lyman-a line. In this paper we thus focus on the 
effects of radiation pressure alone on the structure and the 
spectral signature of the exosphere. 


5.1. Radiation pressure: shaping the cloud 


If the atmosphere of an exoplanet is significantly expanded, 
various regions of the exosphere move on very different 
orbits and revolve around the star in differential rotation. 
In the absence of any radiation pressure from the star, gas 
escaping beyond the Roche Lobe would follow ballistic tra¬ 
jectories subjected to stellar gravity alone. The exosphere 
would suffer shear, the gas within the orbit of the planet 
accelerating and falling toward the star while the gas 
beyond would decelerate and move toward larger distances 
(top panel in Fig. I^. However, this is an unrealistic 
case for neutral hydrogen escape, as radiation pressure 
from the Lyman-a line will always compensate for some 
fraction of the stellar gravity. In the case of close-in planets 
around K-type stars and earlier, such as HD 189733b and 
HD 209458b, the Lyman-a line is actually strong enough 
for radiation pressure to overcome up to 3-5 times the 
gravity, thus swiftly accelerating the H^ atoms away from 
the star at radial velocities exceeding lOOkrns”^ (bottom 
panel in Fig. [^. The exosphere is compressed on its 
dayside, while gas escaping at the limbs of the atmosphere 
is blown away and shaped into a narrow comet ary tail, 
bent toward the radial direction of the stellar photons (e.g., 
Schneider et al.||1998[ |Vidal-Madjar et al.||2QQ3[ [Bourrier 
& Lecavelier des Etangs 2013). Orbiting a M dwarf at 
moderate orbital distancA the warm Neptune GJ 436 
b comparatively receives a lower flux from the stellar 
Lyman-a line. With radiation pressure compensating for 
up to ^60% of the gravitational attraction at low radial 
velocities (Fig. [^, the escaping hydrogen atoms cannot 
be accelerated away from the planet but are nonetheless 
subjected to a strong radiative braking. Indeed, even a 
small reduction of stellar gravity is enough for all the 
exospheric gas to decelerate and move to larger orbits, 
preventing the formation of a stream of matter infalling 
toward the star. In contrast to highly irradiated planets, 
the gas can diffuse within a larger volume around the 
planet and forms a wider cometary tail closer to the planet 
orbit (middle panel in Fig. [^. 


5.2. Effect of radiation pressure on the velocity field 

Radiation pressure has also a major effect on the velocities 
of the escaping hydrogen atoms that correspond to the 
wavelength range of the observed absorption profile. In the 
stellar reference frame, the initial velocity of an escaping 
particle is the combination of the planet orbital velocity 
and the velocity of the planetary wind at the Roche 
Lobe. For close-in exoplanets (a <0.1 an), orbital velocity 
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Fig. 2. Views of GJ436b exosphere (blue dots) from the 
perpendicular to the orbital plane. The planet is the small 
disk at the intersection of the star/Earth LOS (dashed 
black line) and the planetary orbit (green line). All system 
properties are similar in the three panels, except for radia¬ 
tion pressure. Upper panel: Without radiation pressure, 
stellar gravity shears the extended cloud of hydrogen, and 
gas at shorter orbital distance than the planet falls toward 
the star. Middle panel: With even a moderate radiation 
pressure (^60% of stellar gravity at maximum), all the es¬ 
caping gas decelerates and moves to larger orbits, forming 
a longer comet-like tail trailing behind the planet. This is 
the real radiation pressure that corresponds to the observed 
Lyman-(a flux. Lower panel: With a high radiation pres¬ 
sure (^400% of stellar gravity at maximum), the escaping 
gas is swiftly blown away. The planet is surrounded by a 
much smaller coma but trailed by a narrower, more radial, 
cometary-tail. 


is generally higher than 100 km s significantly higher 
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). To a first approximation 


set by the planet orbital velocity Vpi. With even a low 
radiation pressure, hydrogen atoms decelerate with respect 
to the planet. Because radiation pressure does not fully 
overcome the star gravity, hydrogen atoms nevertheless 
remain attracted toward the star and will not move beyond 
the plane tangent to the planet orbit at the time of their 
escape, as illustrated in Fig. As a result of both the 
orbital velocity of the planet and the effect of radiation 
pressure, the projected velocity of the gas thus ranges 
approximately between the projection of Vpi (for atoms 
close to the planet) and -135kms“^, which corresponds 
to the highest planet orbital velocity. Therefore, this 
dynamics is consistent with the detection of blue-shifted 
absorption up to ab out -12Qkms~^ during the transit of 
GJ 436 b exosphere (Ehrenreich et al 


during 1 


5.3. Influence of the physical parameters 

We used simulations performed with EVE to study the 
effects of the most influential star and planet properties 
on the structure of the exosphere and its absorption 
signature. The planetary escape rate and the stellar 
photo-ionization rate are two basic properties of neutral 
hydrogen escape, but in contrast to hot-Jupiters sub¬ 
jected to strong radiation pressure the planetary wind 
velocity also has a significant influence on the exosphere 
of a moderately irradiated, lower-mass planet like GJ 436 b. 


— Planetary wind velocity Upi_wind- This parameter 
influences the spatial and spectral dispersion of the gas 
in the exosphere. Higher values for Upi_wind increase the 
width of the coma and the exospheric tail, extending 
in particular the duration of the pre-transit and post¬ 
transit absorptions. The same amount of gas escaping 
faster from the atmosphere is more diffused and has 
a wider distribution of Doppler velocities, yielding 
lower depths in the core of the absorption profile but 
spreading the absorption over a larger spectral range. 

— Photo-ionization rate Tion- This parameter acts on the 
lifetime of neutral hydrogen atoms, and yields local 
variations in the size of the exosphere and its density. 
Because of self-shielding, photo-ionization has indeed 
more impact on the regions facing the star, which leads 
to a more significant erosion of the tail and the front of 
the exosphere than its denser coma. Higher values for 
Tion thus decrease more strongly the pre-transit and 
post-transit absorption depths. 

— Escape rate M^o: This parameter changes the local 
density roughly uniformly within the exosphere. An 
increase in the escape rate can thus partly balance a 
stronger photo-ionization as long as self-shielding is not 
significant. Higher values for Mho increase the depth of 
the absorption profile at all wavelengths, and lengthen 
the exospheric tail. 
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Fig. 3. Views of GJ 436b exosphere (gray dots) from the perpendicular to the orbital plane, 3.5 h after the optical transit. 
The white dashed line shows the LOS toward Earth. Arrows display the velocity field of the hydrogen atoms with respect 
to the star, colored as a function of the time they escaped the atmosphere (counted from the center of the optical transit). 
In the reduced gravity field from the star, particles move beyond the orbit of the planet. However, radiation pressure does 
not overcome stellar gravity, and particles are still deviated from the tangent to the orbit at the time of their escape (solid 
colored lines). Because radiation pressure is velocity-dependent, the particles do not stay on Keplerian orbits. To better 
illustrate the effect of radiative braking, this simulation was performed with a low ionization rate artificially lengthening 
the exospheric tail. 


When hydrogen densities are high enough, self-shielding 
can protect the interior of the coma from radiation pressure. 
Gas in these regions and escaping from the planet dayside 
is subjected to the nominal stellar gravity, and can move far 
ahead of the planet toward the star. We note however that 
the very high escape rate/very low photo-ionization rate 
that would lead to the formation of this extended leading 
tail are excluded, as it would occult the star too early com¬ 
pared to the observations. 


6. Interpretation of the observations 

6.1. Qualitative description 


The main features of the Ly man-g line observa ti ons o f 
GJ 436 b were described by Ehrenreich et al. (2015). 
They showed that the flux in the entire red wing, and at 
high velocities in the blue wing, is stable over the transit 


duration and between the different epochs. Visit 1 displays 
a puzzling absorption in the red wing just after the optical 
transit, which is not reproduced in subsequent visits. 
Absorption attributed to the exosphere occurs for all visits 
in a similar range of velocities up to about -120kms“^ in 
the blue wing of the line. Overall the light curves measured 
in this range (Eig. present similar variations over time 
for the three visits: 1) the exospheric transit begins about 
2 h before the optical transit; 2) the maximum absorption 
depth is found at about the time of the optical transit; 3) 
although it decreases after the optical transit, absorption 
may still be observed for several hours afterwards. There 
are nevertheless some differences between the visits in 
the way the absorption varies over time. Visit 1 displays 
smooth variations during the atmospheric ingress and 
egress. The decrease in flux is sharper at the ingress of 
the other visits, in particular in Visit 3. Regarding the 
post-transit phases, absorption decreases very quickly in 
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Fig. 4. Views of GJ 436b exosphere (lower panels) and its theoretical absorption profile, as it would be observed close 
to the star (upper panels). We caution that absorption in the range [-40 ; 20]knis“^ cannot be observed from Earth, 
and that STIS LSF would spread this absorption signal over several pixels. Time increases from left to right and is 
indicated on each plot. The dashed black lines limit the LOS toward the stellar surface, while the dotted line indicates 
the star-planet axis. Black arrows show the velocity field of the atoms relative to the star. Gas is seen in the orbital 
plane and colored as a function of its velocity projected on the LOS. Left panel: As they circumvent the planet to move 
away from the star, atoms in the front of the coma absorb the stellar fiux at blueshifted velocities before the optical 
transit. Middle panel: At the time of the optical transit the strongest absorption comes from the core of the coma, and is 
centered on the projection of the planet orbital velocity at -16kms“^. The outer regions of the exospheric tail contribute 
to absorption at higher velocities in the blue wing. Right panel: After the optical transit, absorption in the blue wing of 
the line comes exclusively from the tail. 


Visit 2 while Visit 3 displays absorption depths similar to 
those of the optical transit. 


6.2. Possible scenarios 

In this section we investigate different physical scen arios 
that could explain the features described in Sect. |6.1[ 


6.2.1. Radiative braking 


As shown in Sect. |5.2[ the velocities expected from the 
dynamics of the radiatively braked exosphere are con¬ 
sistent with the range of the absorption signature u 


to -120kms“^. Using the model described in Sect, 
we further found that the combination of low radiation 
pressure with high planetary wind velocities of about 
30kms“^ allows the diffusion of the escaping gas within 
a coma large enough to explain the deep absorption 
observed in all visits. This wind velocity, higher than those 


obtained from the modeling of hot-Jupiters expanding 

thermosphere (^1-10kms see, e.g.,; 

Murray-Clay et al. 

20091 Koskinen et al. 

2013a), may be 

inked to the lower 


a magnetically driven planetary wind (Tanaka et al. 


The extended coma also produces a deeper and earlier 
ingress than observed, in particular for Visit 3 (Fig. [^. 
We also found that the photo-ionization of the gas and its 
natural dilution as it moves away from the planet leads to 
a stronger decrease of the absorption depth than observed 
during the post-transit phases in Visit 3. Therefore, while 
the exosphere is always subjected to radiation pressure and 
this force naturally reproduces the observed velocities of 
the escaping gas, radiative braking alone does not explain 
the time variations of the measured absorption depth. In 
other words, radiation pressure is the main driver of the 
velocity field in the exosphere, and secondary mechanisms 
must reshape the geometry of the gas cloud. 


6.2.2. Radiative blow-out 

As is the case for hot Jupiters, a higher radiation pressure 
could repel the escaping hydrogen and reduce the front 
of the coma responsible for an extended early ingress. 
However, because radiation pressure is proportional to the 
stellar flux, this would require a higher flux in the entire 
Lyman-(a line, implying that the spectra measured during 
the first orbits in Visits 2 and 3 are already absorbed by 
the exosphere. This is unlikely, because these first orbit 
spectra are remarkably similar (as would be the case if 
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Fig. 5. Transit light curves of the Lyman-a line in Visit 1 
(upper panel), Visit 2 (middle panel), and Visit 3 (lower 
panel), integrated between -120 and -40kms“^. Black 
points correspond to the observations, while the black line 
shows the optical transit. The blue lines corresponds to 
best-fit simulations when the exosphere is subjected to stel¬ 
lar gravity and radiation pressure alone. 


they corresponded to the stable intrinsic stellar Lyman-a 
line), in contrast to other spectra measured at similar 
phases of Visits 2 and 3 (Fig. |^. Furthermore to explain 
absorption as early as the first orbits in Visits 2 and 3, 
the coma should be even more extended than is the case 
with a low radiation pressure, in contradiction with its 
narrowing from the stronger radiation pressure. 


6.2.3. Transition region 


Ehrenreich et al. (2015) used the EVE code with similar 
settings as those described in Sect. to fit the three 
combined visits of GJ 436 b in the frame of a radiative 
braking scenario. They reproduced the observed sharp 
early ingress by introducing an additional free parameter 


to the model, which is a multiplicative factor to artificially 
increase the ionization cross-section of atoms. The 
corresponding increase in self-shielding coupled with high 
ionization rates allows for the thin outer envelop of the 
exosphere facing the star to be eroded, while preserving the 
denser inner regions. The sharp density transition region 
thus created is akin to the boundary of a H II region, or 
to a bow-shock environment. However this model requires 
extremely high ionization rates in the order of 10“^ s“^, 
and makes strong assumptions on self-shielding that affect 
the overall structure of the exosphere. Detailed models of 
magnetic interactions between the star and planet will be 
required to investigate this scenario further. 


6.2.4. Stellar wind interactions 


Charge exchange between the escaping atoms and stel¬ 
lar p rotons could also ionize the front layers of the coma 
(e.g., Holmstrom et al. 2008) and reduce the depth of the 
early ingress. These interactions would also ionize planetary 
neutral atoms in every region of the exosphere, decreas¬ 
ing the overall absorption, but this might be compensated 
for by the concomitant creation of an independent popu¬ 
lation of neutralized protons repopulating the exosphere. 
Moving with the high radial velocity of the stellar wind, 
this population would be less sensitive to photo-ionization 
and spatial dilution than the planetary atoms and would 
keep the same distribution over longer time scales, possi¬ 
bly expl ainin g the post-transit absorption depths in Visit 
3 (Sect. |6.1[ ). A combination of the stable radiation pres¬ 
sure with a variable stellar wind may explain the different 
variations in absorption depths observed between the three 
visits. 


7. Discussion 

New HST observations of the Neptune-mass exopla net GJ 
436b revealed a spectacular atmosp heric escape (Kulow 


et al. 20141 Ehrenreich et al. 2015), repeatable over the 


three epochs of observations. Where the exospheres of 
known evaporating hot Jupiters cover ^15% of the stellar 
disk, the cloud of hydrogen surrounding GJ 436 occults up 
to ^60% of the star with significant absorption observed 
several hours before and after the optical transit. The 
amplitude of the absorption and the corresponding size 
of the exosphere may appear surprising, considering that 
GJ 436b is only moderately irradiated by its M-dwarf host 
star. 

To understand this puzzle, we have adapted and run EVE, 
our 3D particle simulation of atmospheric escape, with the 
specific properties of the GJ 436 system. In this paper we 
focused on a scenario where escaping atoms are submitted 
to radiation pressure, stellar and planet gravitational 
forces, and photo-ionization to assess the impact of these 
phenomenon on the observed Lyman-a transit signature. 
Radiation pressure is a velocity-dependent force that can 
be calculated directly from observations of the Lyman-a 
line, and we provide detailed reconstructions of the 
intrinsic stellar line at the different epochs. We found 
that radiation pressure from GJ436 is too low to exceed 
stellar gravity and to repel the escaping hydrogen atoms 
away from the star, as is the case for hot Jupiters. Even 
this moderate radiation pressure is nevertheless enough 
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to brake the gravitational deviation of the gas toward the 
star, and it allows the formation of a large coma comoving 
with the planet (although not gravitationally bound to it) 
that extends in a broad cometary tail trailing the planet 
and progressively eroded by photo-ionization. 

Our simulations show that the low radiation pressure 
from GJ436 explains the velocity range of the observed 
absorption signatures well up to about -120kms“^ and 
suggest that a high planetary wind velocity of about 
30kms“^ is required to create a sizeable coma reproducing 
the observed absorption depths. However, we found that 
the variations in the depth and duration of the absorption 
at the different phases of the transit, some specific to the 
visit, cannot be explained by radiation pressure alone. In 
the low radiation pressure field of the star, the coma also 
extends far ahead of the planet, producing a deeper and 
earlier ingress than observed. An additional mechanism, 
such as interactions with the stellar wind, is thus needed 
to reduce the size of the coma front. By creating a new 
population of neutral hydrogen atoms with the properties 
of the stellar wind, charge exchange may also explain 
peculiar features observed in the last visit. 

We emphasize that future studies intent on reproducing 
the observations of GJ 436 b should take the entire spectral 
information into account and not only the photometric 
light curves. Given the size of the exosphere, the regions 
that were observed transiting the star only sample a small 
part of the atmospheric structure, and more observations 
will be required to fully constrain the dynamics of the 
escaping gas and the properties of GJ 436 b environment. 
HST observations of a partially tr ansiting exosphere 
around the warm giant 55 Cancri b (jEhrenreich et al. 


2012) previously hinted that moderately irradiated planets 


could have significantly extended exospheres. The detec¬ 
tion of GJ 436b exosphere offers the opportunity to probe 
the atmospheric escape of a “cool” exoplanet in detail 
(^800 K), whose atmospheric photochemistry m ay also 


be significa ntly affected by Lyman-a radiation (Miguel 
|et al. [ 2 M 5 ). Generalizing these results to other systems, 
in particular those that will be detected around M dwarfs 
by CHEOPS, TESS, and PLATO, will open tantalizing 
new perspectives for characterizating warm and temperate 
exoplanets. 
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